The R Spatial Ecosystem

A few main packages are dedicated to spatial data handling like import/export, display, geoprocessing, spatial ststistics.

“Historical” Packages

rgdal: interface between R and GDAL (Geospatial Data Abstraction Library) and PROJ4 libraries.
sp: classes and methods for spatial data in R.
rgeos: interface between GEOS (Geometry Engine - Open Source) library and R: area, perimeter, distances, dissolve, buffer, overlap, union, contains…

These packages are still widely used.

Simple Features for R

sp, rgeos and rgdal functionnalities in one package.
Easier data handling, simpler objects.
Main author and maintainer: Edzer Pebesma (also sp author)
Compatibility with the pipe synthax and tidyverse operators.
First release: October 20, 2016
Website: Simple Features for R

sf objects data structure:

format sf

Using `sf`

Data Import

library(sf)
mtq <- st_read("../data/mtq/martinique.shp")

Reading layer `martinique' from data source `/home/tim/Documents/prz/flo/data/mtq/martinique.shp' using driver `ESRI Shapefile'
Simple feature collection with 34 features and 23 fields
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 690574.4 ymin: 1592426 xmax: 736126.5 ymax: 1645660
epsg (SRID):    32620
proj4string:    +proj=utm +zone=20 +datum=WGS84 +units=m +no_defs

Data Display

Default

plot(mtq)

Only geometry

plot(st_geometry(mtq))

Centroids Extraction

mtq_c <- st_centroid(mtq)
plot(st_geometry(mtq))
plot(st_geometry(mtq_c), add = TRUE, cex = 1.2, col = "red", pch = 20)

Distance Matrix

mat <- st_distance(x = mtq_c, y = mtq_c)
mat[1:5, 1:5]

Units: m
          [,1]     [,2]      [,3]      [,4]      [,5]
[1,]     0.000 35297.56  3091.501 12131.617 17136.310
[2,] 35297.557     0.00 38332.602 25518.913 18605.249
[3,]  3091.501 38332.60     0.000 15094.702 20226.198
[4,] 12131.617 25518.91 15094.702     0.000  7177.011
[5,] 17136.310 18605.25 20226.198  7177.011     0.000

Polygons Aggregation

mtq_u <- st_union(mtq)
plot(st_geometry(mtq), col = "lightblue")
plot(st_geometry(mtq_u), add = T, lwd = 2, border = "red")

Buffer Zone

mtq_b <- st_buffer(x = mtq_u, dist = 5000)
plot(st_geometry(mtq), col = "lightblue")
plot(st_geometry(mtq_u), add = T, lwd = 2)
plot(st_geometry(mtq_b), add = T, lwd = 2, border = "red")

Polygon Intersection

m <- rbind(c(700015, 1624212), c(700015, 1641586), c(719127, 1641586), c(719127, 
    1624212), c(700015, 1624212))
p <- st_sf(st_sfc(st_polygon(list(m))), crs = st_crs(mtq))
plot(st_geometry(mtq))
plot(p, border = "red", lwd = 2, add = T)

mtq_z <- st_intersection(x = mtq, y = p)
plot(st_geometry(mtq))
plot(st_geometry(mtq_z), col = "red", border = "green", add = T)

Voronoi Polygons

google: “st_voronoi R sf” (https://github.com/r-spatial/sf/issues/474 & https://stackoverflow.com/questions/45719790/create-voronoi-polygon-with-simple-feature-in-r)

mtq_v <- st_voronoi(x = st_union(mtq_c))
mtq_v <- st_intersection(st_cast(mtq_v), st_union(mtq))
mtq_v <- st_join(x = st_sf(mtq_v), y = mtq_c, join = st_intersects)
mtq_v <- st_cast(mtq_v, "MULTIPOLYGON")
plot(st_geometry(mtq_v), col = "lightblue")

Maps with R

Several solutions are available:

ggplot2 users can have a look to ggplot2 mapping features (geom_sf) that can mix nicely with ggspatial.
For more advanced mapping features in a ggplot2-like syntax have a look to tmap
cartography is based on base graphics and allow most of basic and advanced cartographic representations. Full disclosure: the speaker is the maintainer of cartography.

Here we will focus on cartography and do small examples with ggplot2.

`cartography`

Data Preparation

library(sf)
# Import geo layer
sm <- st_read(dsn = "../data/rhone/seine_maritime.geojson", stringsAsFactors = F)

Reading layer `OGRGeoJSON' from data source `/home/tim/Documents/prz/flo/data/rhone/seine_maritime.geojson' using driver `GeoJSON'
Simple feature collection with 718 features and 1 field
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: 0.06560943 ymin: 49.25139 xmax: 1.790225 ymax: 50.07121
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs

dep <- st_read(dsn = "../data/rhone/dep.geojson", stringsAsFactors = F)

Reading layer `OGRGeoJSON' from data source `/home/tim/Documents/prz/flo/data/rhone/dep.geojson' using driver `GeoJSON'
Simple feature collection with 96 features and 2 fields
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: -5.138996 ymin: 41.36216 xmax: 9.559598 ymax: 51.089
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs

# change projection (lambert93)
sm <- st_transform(sm, 2154)
dep <- st_transform(dep, 2154)

# Import dataset
csp <- read.csv("../data/rhone/data_seine_maritime.csv")
# merge geolayer and dataset
sm <- merge(sm, csp, by = "INSEE_COM", all.x = TRUE)

Proportional symbols

library(cartography)

plot(st_geometry(sm))
propSymbolsLayer(sm, var = "act")
title("Active Population")

# Custom map of active population
par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
bb <- st_bbox(sm)
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
plot(st_geometry(sm), col = "cornsilk2", border = NA, lwd = 0.5, add = T)
propSymbolsLayer(sm, var = "act", col = "darkblue", inches = 0.5, border = "white", 
    lwd = 0.7, symbols = "square", legend.style = "e", legend.pos = "topleft", 
    legend.title.txt = "Labor Force\n(2014)", legend.values.rnd = 0)
barscale(size = 10)
north(pos = "topright", col = "darkblue")
layoutLayer(title = "Workforce in Seine-Maritime", sources = "Insee, 2018", 
    author = "Kim & Tim, 2018", col = "darkblue", coltitle = "white", frame = TRUE, 
    scale = NULL, north = FALSE)

Typology Map

par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
typoLayer(sm, var = "cat", border = "ivory", lwd = 0.5, legend.values.order = c("agr", 
    "art", "cad", "int", "emp", "ouv"), col = c("#e3b4a2", "#a2d5d6", "#debbd4", 
    "#b5dab6", "#afc2e3", "#e9e2c1"), add = TRUE, legend.pos = "n")

legendTypo(title.txt = "Dominant Category", col = c("#e3b4a2", "#a2d5d6", "#debbd4", 
    "#b5dab6", "#afc2e3", "#e9e2c1"), categ = c("Agriculteurs", "Artisans", 
    "Cadres", "Prof. Inter.", "Employés", "Ouvriers"), nodata = F)
barscale(size = 10)
north(pos = "topright", col = "darkblue")
layoutLayer(title = "Workforce Distribution in Seine-Maritime", sources = "Insee, 2018", 
    author = "Kim & Tim, 2018", col = "darkblue", coltitle = "white", frame = TRUE, 
    scale = NULL, north = FALSE)

Choropleth Map

sm$pcad <- 100 * sm$cad/sm$act

bks <- getBreaks(v = sm$pcad, method = "quantile", nclass = 6)
cols <- carto.pal("green.pal", 3, "wine.pal", 3)
par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
choroLayer(sm, var = "pcad", breaks = bks, col = cols, border = "grey80", lwd = 0.4, 
    legend.pos = "topleft", legend.title.txt = "Share of managers\n(%)", add = T)
layoutLayer(title = "Managers", sources = "Insee, 2018", author = "Kim & Tim, 2018", 
    theme = "green.pal", col = "darkred", coltitle = "white", postitle = "center", 
    frame = TRUE, scale = 10)
north(pos = "topright", south = TRUE)

library(cartogram)
sm_c1 <- cartogram_cont(sm, "act", 3)
sm$pcad <- 100 * sm$cad/sm$act

bks <- getBreaks(v = sm$pcad, method = "quantile", nclass = 6)
cols <- carto.pal("green.pal", 3, "wine.pal", 3)
par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
choroLayer(sm_c1, var = "pcad", breaks = bks, col = cols, border = "grey80", 
    lwd = 0.4, legend.pos = "topleft", legend.title.txt = "Share of managers\n(en %)", 
    add = T)
layoutLayer(title = "Managers", sources = "Insee, 2018", author = "Kim & Tim, 2018", 
    theme = "green.pal", col = "darkred", coltitle = "white", postitle = "center", 
    frame = TRUE, scale = 10)
north(pos = "topright", south = TRUE)

Gridded Map

These maps may allow to hide or, at least, diminish the MAUP.

sqrt(quantile(st_area(sm), 0.9))

3904.689 m

grid <- getGridLayer(x = sm, cellsize = 3900 * 3900, type = "regular", var = c("cad", 
    "act"))
grid$pcad <- 100 * grid$cad/grid$act
par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
choroLayer(grid, var = "pcad", breaks = bks, col = cols, border = "grey80", 
    lwd = 0.4, legend.pos = "topleft", legend.title.txt = "Share of managers\n(en %)", 
    add = T)
layoutLayer(title = "Managers", sources = "Insee, 2018", author = "Kim & Tim, 2018", 
    theme = "green.pal", col = "darkred", coltitle = "white", postitle = "center", 
    frame = TRUE, scale = 10)
north(pos = "topright", south = TRUE)

Smoothed Map

grid$cad100 <- grid$cad * 100
par(mar = c(0.2, 0.2, 1.4, 0.2), bg = "azure")
plot(st_geometry(dep), col = "ivory", border = "ivory3", xlim = bb[c(1, 3)], 
    ylim = bb[c(2, 4)])
smoothLayer(x = grid, var = "cad100", var2 = "act", typefct = "exponential", 
    span = 4000, beta = 2, breaks = bks, col = cols, legend.pos = "topleft", 
    mask = st_buffer(sm, 0), legend.title.txt = "Share of managers\n(%)*", border = "grey90", 
    lwd = 0.2, add = T)
layoutLayer(title = "Managers", sources = "Insee, 2018", author = "Kim & Tim, 2018", 
    theme = "green.pal", col = "darkred", coltitle = "white", postitle = "center", 
    frame = TRUE, scale = 10)
north(pos = "topright", south = TRUE)
text(x = 518000, y = 6921123, font = 3, cex = 0.8, labels = "*Lissage par potentiels\n fonction exponnentielle\n span = 2.5km, beta = 2")

`ggplot2`

Proportional symbols layer:

library(ggplot2)

ggplot() +
  geom_sf(data = dep, colour = "ivory3",fill = "ivory") +
  geom_sf(data = sm %>%  st_centroid(),
          
          aes(size= cad), colour="#E84923CC", show.legend = 'point') +
  scale_size(name = "Active population",
             breaks = c(100,1000,5000),
             range = c(0,20)) +
  coord_sf(crs = 2154, datum = NA,
           xlim = st_bbox(sm)[c(1,3)],
           ylim = st_bbox(sm)[c(2,4)]
  ) +
  theme_minimal() +
  theme(panel.background = element_rect(fill = "azure",color=NA)) +
  labs(title = "Active population",
       caption = "Insee, 2018\nKim & Tim, 2018")

Choropleth layer:

ggplot() +
  geom_sf(data = dep, colour = "ivory3",fill = "ivory") +
  geom_sf(data = sm, aes(fill = pcad), colour = "grey80") +
  scale_fill_gradientn(name = "Share of managers\n(%)",
                       
                       colours = carto.pal("green.pal", 3,"wine.pal",3),
                       values=bks/max(bks))+
  coord_sf(crs = 2154, datum = NA,
           xlim = st_bbox(sm)[c(1,3)],
           ylim = st_bbox(sm)[c(2,4)]
  ) +
  theme_minimal() +
  theme(panel.background = element_rect(fill = "azure",color=NA)) +
  labs(title = "Managers",
       caption = "Insee, 2018\nKim & Tim, 2018")

Exploration of sirene & osm db for restaurants

CRAN Task View: Analysis of Spatial Data

See file data_prep.R for data extraction.

# spatial data management
library(sf)
# cartography
library(cartography)
# smooth density computation
library(spatstat)
library(raster)
library(maptools)
# interactive mapping
library(mapview)

Display Points

# load data
adm <- readRDS("../data/iris_p13.rds")
feat_sir <- readRDS("../data/sir_p13.rds")
feat_osm <- readRDS("../data/osm_p13.rds")

# Define margins
par(mar = c(0, 0, 1.2, 0), mfrow = c(1, 2))

# plot the maps
plot(st_geometry(adm), col = "ivory1", border = "ivory3", lwd = 0.5, bg = "#FBEDDA")
plot(st_geometry(feat_sir), col = "red", pch = 20, add = T, cex = 0.5)
layoutLayer(title = "SIR", sources = paste0(nrow(feat_sir), " restaurants"), 
    author = "Kim & Tim, 2018", scale = NULL, frame = FALSE, tabtitle = TRUE)
legend(x = "topright", legend = "= 1 restaurant", pch = 20, pt.cex = 0.5, bty = "n", 
    cex = 0.7, col = "red")

plot(st_geometry(adm), col = "ivory1", border = "ivory3", lwd = 0.5, bg = "#FBEDDA")
plot(st_geometry(feat_osm), col = "red", pch = 20, add = T, cex = 0.5)
layoutLayer(title = "OSM", sources = paste0(nrow(feat_osm), " restaurants"), 
    author = "Kim & Tim, 2018", scale = 0.5, north = T, frame = FALSE, tabtitle = TRUE)

Interactive Visualisation

Explore datasets interactively

The default is not bad:

mapview(feat_sir) + mapview(feat_osm, col.regions = "red")

mapview allow hightly customizable maps:

mapview(feat_sir, map.types = "OpenStreetMap", col.regions = "#940000", label = paste(feat_sir$L2_NORMALISEE, 
    feat_sir$NOMEN_LONG, sep = " - "), color = "white", legend = TRUE, layer.name = "SIRENE", 
    homebutton = FALSE, lwd = 0.5) + mapview(feat_osm, col.regions = "#000094", 
    color = "white", legend = TRUE, label = feat_osm$name, lwd = 0.5, layer.name = "OSM", 
    homebutton = FALSE)

These maps, as appealing as they seem to be, are not really suitable for displaying geostatistical info. Nonethless they can be really useful for exploratory data analisys.

Count Points in Units

In admin units

# intersect adm and feat
inter_osm <- st_intersects(adm, feat_osm)
inter_sir <- st_intersects(adm, feat_sir)
# count points in polygons
adm <- st_sf(adm[, 1:2, drop = T], n_osm = sapply(X = inter_osm, FUN = length), 
    n_sir = sapply(X = inter_sir, FUN = length), geometry = st_geometry(adm))
# Display the map
par(mar = c(0, 0, 1.2, 0), mfrow = c(1, 2))
plot(st_geometry(adm), col = "ivory1", border = "ivory3", lwd = 0.5, bg = "#FBEDDA")
propSymbolsLayer(adm, var = "n_sir", inches = 0.1)
layoutLayer(title = "SIR", sources = "", author = "", scale = NULL, tabtitle = TRUE, 
    frame = FALSE)
# Display the map
plot(st_geometry(adm), col = "ivory1", border = "ivory3", lwd = 0.5, bg = "#FBEDDA")
propSymbolsLayer(adm, var = "n_osm", inches = 0.1)
layoutLayer(title = "OSM", sources = "", author = "", scale = 0.5, north = T, 
    tabtitle = TRUE, frame = FALSE)

In a regular grid

Create a small function that build a regular grid and count points in each cells.

pt_in_grid <- function(feat, adm, cellsize = 50) {
    # create a regular grid (adm bbox)
    grid <- st_make_grid(x = adm, cellsize = cellsize, what = "polygons")
    
    # keep only cells that intersects adm
    . <- st_intersects(grid, adm)
    grid <- grid[sapply(X = ., FUN = length) > 0]
    
    # count pts in grid
    . <- st_intersects(grid, feat)
    grid <- st_sf(n = sapply(X = ., FUN = length), grid)
    
    return(grid)
}

Create a function that plot a map

plot_grid <- function(grid, adm, title, bks, col) {
    plot(st_geometry(adm), col = "ivory1", border = "ivory3", lwd = 0.5)
    choroLayer(grid, var = "n", border = NA, add = TRUE, breaks = bks, col = cols, 
        legend.pos = "n")
    layoutLayer(title = title, scale = NULL, frame = FALSE, tabtitle = TRUE, 
        author = "", sources = "")
}

Builds grids and display maps

# plot 6 maps on a single figure
par(mar = c(0, 0, 1.8, 0), mfrow = c(2, 3), bg = "#FBEDDA", ps = 16)

# defines a unique set of breaks for all the maps
bks <- seq(1, 22, 3)
cols <- carto.pal("sand.pal", length(bks) - 1)

# SIR maps 50 m cells
grid <- pt_in_grid(feat_sir, adm, 50)
plot_grid(grid = grid, adm = adm, title = "SIR - 50 m cells", col = cols, bks = bks)

# 100 m cells
grid <- pt_in_grid(feat_sir, adm, 100)
plot_grid(grid = grid, adm = adm, title = "SIR - 100 m cells", col = cols, bks = bks)

# 150 m cells
grid <- pt_in_grid(feat_sir, adm, 150)
plot_grid(grid = grid, adm = adm, title = "SIR - 150 m cells", col = cols, bks = bks)

# Add a north arrow
north()

# OSM maps 50 m cells
grid <- pt_in_grid(feat_osm, adm, 50)
plot_grid(grid = grid, adm = adm, title = "OSM - 50 m cells", col = cols, bks = bks)

# Add sources
mtext(text = "Map data © OpenStreetMap contributors, under CC BY SA.\nSIRENE - 2018, INSEE - 2018", 
    side = 1, line = -1, cex = 0.5, adj = 0)

# 100 m cells
grid <- pt_in_grid(feat_osm, adm, 100)
plot_grid(grid = grid, adm = adm, title = "OSM - 100 m cells", col = cols, bks = bks)

# 150 cells
grid <- pt_in_grid(feat_osm, adm, 150)
plot_grid(grid = grid, adm = adm, title = "OSM - 150 m cells", col = cols, bks = bks)

# add legend
legendChoro(pos = "topright", title.cex = 0.7, values.cex = 0.5, title.txt = "N. restaurants", 
    breaks = bks, nodata = FALSE, values.rnd = 0, col = cols)

# add scalebar
barscale(1)

Smooth Density Analysis

Define a function that compute smoothed restaurant density (KDE method).

compute_kde <- function(feat, adm, title, sigma = 100, res = 50) {
    # define an observation window
    w <- as.owin(as(adm, "Spatial"))
    # sf to coords
    pts <- st_coordinates(feat)
    # coord to ppp
    p <- ppp(pts[, 1], pts[, 2], window = w)
    # Compute KDE
    dens <- density.ppp(p, sigma = sigma, eps = res)
    # im to raster (+ proj & km2)
    result <- raster(dens, crs = st_crs(adm)[[2]]) * 1e+06
    return(result)
}

Define a function that map the results

plot_kde <- function(x, adm, method = "equal", nclass = 12, title = "", author = "", 
    sources = "", scale = NULL, legend.title = "Density\n(n/km2)", legend.pos = "topright", 
    values.rnd = 0) {
    # compute breaks
    bks <- unique(getBreaks(values(x), nclass = nclass, method = method))
    # color ramp
    cols <- (mapview::mapviewGetOption("raster.palette"))(length(bks) - 1)
    # plot the map
    plot(st_geometry(adm), col = NA, border = NA, main = "", bg = "#FBEDDA")
    plot(x, breaks = bks, col = cols, add = T, legend = F)
    legendChoro(pos = legend.pos, title.txt = legend.title, values.cex = 0.5, 
        breaks = bks, nodata = FALSE, values.rnd = values.rnd, col = cols)
    layoutLayer(title = title, scale = scale, tabtitle = TRUE, frame = FALSE, 
        author = author, sources = sources)
}

# plot 6 maps on a single figure
par(mar = c(1, 0, 1.8, 0), mfrow = c(2, 3), bg = "#FBEDDA", ps = 16)

x <- compute_kde(feat = feat_sir, adm = adm, sigma = 50, res = 20)
plot_kde(x, adm, title = "SIR - Sigma = 50", legend.title = "Restaurants\n/km2")
x <- compute_kde(feat = feat_sir, adm = adm, sigma = 100, res = 20)
plot_kde(x, adm, title = "SIR - Sigma = 100", legend.title = "Restaurants\n/km2")
x <- compute_kde(feat = feat_sir, adm = adm, sigma = 250, res = 20)
plot_kde(x, adm, title = "SIR - Sigma = 250", legend.title = "Restaurants\n/km2")
x <- compute_kde(feat = feat_osm, adm = adm, sigma = 50, res = 20)
plot_kde(x, adm, title = "OSM - Sigma = 50", legend.title = "Restaurants\n/km2")
x <- compute_kde(feat = feat_osm, adm = adm, sigma = 100, res = 20)
plot_kde(x, adm, title = "OSM - Sigma = 100", legend.title = "Restaurants\n/km2")
x <- compute_kde(feat = feat_osm, adm = adm, sigma = 250, res = 20)
plot_kde(x, adm, title = "SIR - Sigma = 250", legend.title = "Restaurants\n/km2")

reproducibility

sessionInfo()

R version 3.5.1 (2018-07-02)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Debian GNU/Linux 9 (stretch)

Matrix products: default
BLAS: /usr/lib/libblas/libblas.so.3.7.0
LAPACK: /usr/lib/lapack/liblapack.so.3.7.0

locale:
 [1] LC_CTYPE=fr_FR.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=fr_FR.UTF-8        LC_COLLATE=fr_FR.UTF-8    
 [5] LC_MONETARY=fr_FR.UTF-8    LC_MESSAGES=fr_FR.UTF-8   
 [7] LC_PAPER=fr_FR.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] mapview_2.5.0       maptools_0.9-3      raster_2.6-7       
 [4] sp_1.3-1            spatstat_1.56-1     rpart_4.1-13       
 [7] nlme_3.1-137        spatstat.data_1.3-1 ggplot2_3.0.0      
[10] cartography_2.1.2   sf_0.6-3            rmdformats_0.3.3   
[13] knitr_1.20         

loaded via a namespace (and not attached):
 [1] viridisLite_0.3.0    shiny_1.1.0          assertthat_0.2.0    
 [4] highr_0.7            stats4_3.5.1         yaml_2.2.0          
 [7] pillar_1.3.0         backports_1.1.2      lattice_0.20-35     
[10] glue_1.3.0           digest_0.6.16        promises_1.0.1      
[13] polyclip_1.9-1       colorspace_1.3-2     htmltools_0.3.6     
[16] httpuv_1.4.5         Matrix_1.2-14        plyr_1.8.4          
[19] pkgconfig_2.0.2      questionr_0.6.3      bookdown_0.7        
[22] purrr_0.2.5          xtable_1.8-2         scales_1.0.0        
[25] webshot_0.5.0        tensor_1.5           satellite_1.0.1     
[28] later_0.7.3          spatstat.utils_1.9-0 tibble_1.4.2        
[31] mgcv_1.8-24          withr_2.1.2          lazyeval_0.2.1      
[34] magrittr_1.5         crayon_1.3.4         mime_0.5            
[37] deldir_0.1-15        evaluate_0.11        foreign_0.8-70      
[40] class_7.3-14         tools_3.5.1          formatR_1.5         
[43] stringr_1.3.1        munsell_0.5.0        bindrcpp_0.2.2      
[46] compiler_3.5.1       e1071_1.7-0          rlang_0.2.2         
[49] classInt_0.2-3       units_0.6-0          grid_3.5.1          
[52] rstudioapi_0.7       htmlwidgets_1.2      goftest_1.1-1       
[55] crosstalk_1.0.0      miniUI_0.1.1.1       base64enc_0.1-3     
[58] labeling_0.3         rmarkdown_1.10       gtable_0.2.0        
[61] abind_1.4-5          DBI_1.0.0            R6_2.2.2            
[64] dplyr_0.7.6          rgeos_0.3-28         bindr_0.1.1         
[67] rprojroot_1.3-2      stringi_1.2.4        Rcpp_0.12.18        
[70] png_0.1-7            leaflet_2.0.1        spData_0.2.9.3      
[73] tidyselect_0.2.4     xfun_0.3

Geovisualization and Spatial Analysis